# aem.single.line_sink --
#
# Copyright (c) 2000, V.Kelson and M.Bakker

from types import *
from cmath import log,pi

import geometry
XY = geometry.xy.XY
XYZ = geometry.xyz.XYZ

import base

import aem.container
RowList = aem.container.RowList
ColumnList = aem.container.ColumnList
TraceCheck = aem.container.TraceCheck

class LineSink(geometry.two_d.Line,base.Element):
	'''Class for given constant strength linesinks

	Python call: LineSink(pt1,pt2,sigma,parent)

	Attributes:

	  - pt1: first XY of line (from Line)
	  - pt2: second XY of line (from Line)
	  - width: width of line (from Line)
	  - z1: complex coordinate of pt1 (from Line)
	  - z2: complex coordinate of pt2 (from Line)
	  - zc: complex coordinate of center of line (from Line)
	  - zl: complex coordinate 0.5(z2-z1) (from Line)
	  - halflength: half the length of the line (from Line)
	  - parent: parent of the element (inherited from Element)
	  - parameters: ColumnList of strength parameters, set to sigma (inherited from Element)
	  - solver: solver (inherited from Element)
	  - elementList: list of elements, defaults to itself (inherited from Element)
	  - __numberOfParameters: number of parameters, set to 1
	  - factor: length/4pi
	  - logzl: log(zl)

	Overloaded functions from Element: __init__, potentialInfluence, dischargeInfluence
	'''

	def __init__(self,pt1=XY(-1.0,0.0),pt2=(1.0,0.0),sigma=0.0,parent=None):
		'''Constructor'''
		# Call base class constructors
		geometry.two_d.Line.__init__(self,pt1,pt2)
		base.Element.__init__(self,parent)
		self.__numberOfParameters = 1
		self.parameters = ColumnList([ sigma ])
		self.factor = self.length()/(4.0*pi)
		self.logzl = log(self.zl)
		self.width = 0.1

	def __repr__(self):
		# Parent not yet included
		return 'LineSink' + str((self.pt1,self.pt2,self.parameters[0]))

	def __complexPotentialInfluence(self,xy):
		Z = self.map(xy)
		tol = self.width/self.halfLength
		if (abs(Z+1.0)<tol):
			Z=Z-tol
		elif (abs(Z-1.0)<tol):
			Z=Z+tol
		om1 = (Z+1.0)*log(Z+1.0)
		om2 = (Z-1.0)*log(Z-1.0)
		return self.factor * (om1 - om2 + 2.0*self.logzl - 2.0)

	def potentialInfluence(self,xy,t=0):
		return RowList([(self.__complexPotentialInfluence(xy)).real])

	def __complexDischargeInfluence(self,xy):
		Z = self.map(xy)
		return -self.factor*log((Z+1)/(Z-1))/self.zl

	def dischargeInfluence(self,xy,t=0):
		comdis = self.__complexDischargeInfluence(xy)
		return RowList([XY(comdis.real,-comdis.imag)])

	def nearElement(self,xyz1,xyz2,step):
		xy1 = xyz1.XY(); Z1 = self.map(xy1)
		xy2 = xyz2.XY(); Z2 = self.map(xy2)
		rv = TraceCheck()
		# If point 1 on one side and point 2 on other side
		if (Z1.imag > 0 and Z2.imag < 0) or (Z1.imag < 0 and Z2.imag > 0):
			Xintersec = Z1.real + (0-Z1.imag)/(Z2.imag-Z1.imag) * (Z2.real-Z1.real)
			if abs(Xintersec) <= 1:
				rv.change = 1
				if self.parameters[0] > 0: # Should be modified when normalDischarge is implemented
					xynew = self.mapback(complex(Xintersec,0))
					znew = xyz1.z + (0-Z1.imag)/(Z2.imag-Z1.imag) * (xyz2.z-xyz1.z)
					xyzchanged = XYZ(xynew.x,xynew.y,znew)
					rv.stop = 1
				else: # Should be modified when normalDischarge is implemented
					if Z1.imag < 0:
						xynew = self.mapback(complex(Xintersec,0.001)) # Put 0.001*0.5*L off line-sink
					else:
						xynew = self.mapback(complex(Xintersec,-0.001))
					znew = xyz1.z + (0-Z1.imag)/(Z2.imag-Z1.imag) * (xyz2.z-xyz1.z)
					xyzchanged = XYZ(xynew.x,xynew.y,znew)
				rv.xyzchanged = xyzchanged
		return rv
